package RiskAlert;

import java.util.Arrays;

import RiskAlert.InverseMatrix;

public class GM {

	private static double[] input;//原始数据
	private static double[] GenerateData;//生成数据，累加生成
	private static double a;//发展系数,反映预测的发展态势
	private static double u;//灰作用量,反映数据的发展关系
	private static double[] param=new double[2];//参数矩阵，存放由公式计算得出的参数—a和u，[a,u]


	public GM(double[] input){
		this.input=input;
		GenerateData=new double[input.length];
	}

	/**
	 * 原始数据转换为生成数据
	 * 利用累加方法
	 * 注释掉的部分为计算优化之前的累加方法
	 * @param input
	 * @return
	 */
	private boolean generateData(double[] input){
//		for(int i=0;i<input.length;i++)
//			for(int j=0;j<=i;j++)
//				GenerateData[i]+=input[j];
		GenerateData[0]=input[0];
		for(int i=1;i<input.length;i++)
			GenerateData[i]=GenerateData[i-1]+input[i];
		return true;
	}
	
	/**
	 * 利用矩阵相成求解未知系数a和u
	 * 中国商业银行贷款风险分析.pdf上面的B矩阵有误。参照灰色预测.pdf
	 * @param input
	 * @param GenerateData
	 * @return
	 */
	private boolean solveParam(double[] input,double[] GenerateData){
		int length=input.length-1;
		double[][] B=new double[length][2];
		double[][] Bt=new double[2][length];//B的转置矩阵
		double[][] BtB=new double[length][length];//Bt*B
		double[][] InverseBtB=new double[2][2];//BtB的逆矩阵
		double[][] InverseBtBBt=new double[2][length];//InverseBtB*Bt
		double[][] Y=new double[length][1];
		double[] Yt=Arrays.copyOfRange(input, 1, length);//Y的转置矩阵
		double[][] paramtemp=new double[2][1];
		//初始化Y
		for(int i=0;i<length;i++)
			Y[i][0]=Yt[i];
		//初始化B
		for(int i=0;i<length;i++){
			B[i][0]=-(GenerateData[i]+GenerateData[i+1])/2;
			B[i][1]=1;
		}
		
		/*这个比较繁琐，后面的方法将矩阵乘法写为一个方法，直接调用

		//计算矩阵乘法，Bt*B，真恶心……
		for(int i=0;i<2;i++)//Bt矩阵按行遍历
			for(int j=0;j<2;j++)//B矩阵按列遍历
				for(int k=0;k<length;k++)//Bt的行与B列的乘积
					BtB[i][j]+=Bt[i][k]*B[k][j];
		InverseBtB=InverseMatrix.getIM(BtB);
		//矩阵乘法，InverseBtB*Bt
		for(int i=0;i<2;i++)
			for(int j=0;j<length;j++)
				for(int k=0;k<2;k++)
					InverseBtBBt[i][j]+=InverseBtB[i][k]*Bt[k][j];
		//最后一次InverseBtBBt*Y
		for(int i=0;i<2;i++)
			for(int j=0;j<1;j++)
				for(int k=0;k<length;k++)
					param[i][j]+=InverseBtBBt[i][j]*Y[k][j];
		
		*/
		
		BtB=multiMatrix(Bt,B,2,length,2);//计算BtB
		InverseBtB=InverseMatrix.getIM(BtB);//得到InverseBtB
		
		InverseBtBBt=multiMatrix(InverseBtB,Bt,2,2,length);//计算InverseBtBBt
		
		paramtemp=multiMatrix(InverseBtBBt,Y,2,length,1);

		for(int i=0;i<paramtemp.length;i++)
			param[i]=paramtemp[i][0];
		
		return true;
		
	}
	
	/**
	 * 矩阵乘法
	 * @param n1 第一个矩阵的行数
	 * @param n2 第一个矩阵的列数（第二个矩阵的行数）
	 * @param n3 第二个矩阵的列数
	 * @param matrix1 第一个矩阵，被乘矩阵
	 * @param matrix2 第二个矩阵，乘矩阵
	 * @return 返回两个矩阵相乘结果
	 */
	private double[][] multiMatrix(double[][] matrix1,double[][] matrix2,int n1,int n2,int n3){
		double[][] result=new double[n1][n3];
		for(int i=0;i<n1;i++)
			for(int j=0;j<n3;j++)
				for(int k=0;k<n2;k++)
					result[i][j]+=matrix1[i][k]*matrix2[k][j];
		return result;
	}
	
	/**
	 * 求已经还原的预测值
	 * 中国商业银行贷款风险.pdf（P90）与灰色预测.pdf（P5）在这个公式上面有分歧。这里按照前者公式计算
	 * @param param 参数a，u
	 * @param input 原始数据矩阵的第一个值
	 * @param k 预测值参数
	 * @return 返回第k+1个值的预测值
	 */
	public double getPredictedValue(double[] param,double input,int k){
		//这里有Math.expm1(double x)和Math.exp(double x)两个公式可以选择,这里选用前者，如果对计算影响较大，可以用后者试验一下
		return Math.expm1(param[0])*(input-param[1]/param[0])*Math.exp(-param[0]*k);
	}
	
}
